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ABSTRACT 

Nonlinear kinetic theory of cosmic ray (CR) acceleration in supernova remnants is employed to 
calculate CR spectra. The magnetic field in SNRs is assumed to be significantly amplified by the 
efficiently accelerating nuclear CR component. It is shown that the calculated CR spectra agree in a 
satisfactory way with the existing measurements up to the energy 10 17 eV. The power law spectrum 
of protons extends up to the energy 3 x 10 15 eV with a subsequent exponential cutoff. It gives a 
natural explanation for the observed knee in the Galactic CR spectrum. The maximum energy of the 
accelerated nuclei is proportional to their charge number Z. Therefore the break in the Galactic CR 
spectrum is the result of the contribution of progressively heavier species in the overall CR spectrum so 
that at 10 17 eV the CR spectrum is dominated by iron group nuclei. It is shown that this component 
plus a suitably chosen extragalactic CR component can give a consistent description for the entire 
Galactic CR spectrum. 

Subject headings: cosmic rays - shock waves - acceleration of particles - supernova remnants 
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1. INTRODUCTION 

The main reason why supernova remnants (SNRs) are 
usually considered as the cosmic ray (CR) source is a sim- 
ple argument about the energy required to sustain the 
Galactic cosmic ray (GCR) population against loss by 
escape, nuclear interactions and ionization energy loss. 
Supernovae have enough power to drive the GCR accel- 
eration if there exists a mechanism for channeling about 
10% of the mechanical energy release into relativistic par- 
ticles. The high velocity ejecta produced in the super- 
nova (SN) explosion interacts with the ambient medium 
to produce a strong blast wave. This outer shock may 
accelerate a small suprathermal fraction of the ambient 
plasma to high energies. 

The only theory of particle acceleration which at 
present is sufficiently well developed and specific to allow 
quantitative model calculations is diffusive acceleration 
applied t o the strong out er shock associated with SNRs 
(e.g. see Berczhko 2005, for review). Considerable ef- 
forts have been made during the last years to empirically 
confirm the theoretical expectation that the main part of 
GCRs indeed originates in SNRs. Theoretically, progress 
in the solution of this problem has been due to the devel- 
opment of a kinetic nonlinear theory of diffusive shock ac- 
celera tion (jBerezhko et alJll996t IBerezhko k Volklll997t 
2000) . This theory attempts to include all the most rele- 
vant physical factors, essential for the evolution and CR 
acceleration in a SNR, at least in its early very energetic 
stages, and it is able to make quantitative predictions of 
the expected properties of CRs produced in young SNRs 
and their nonthermal radiation. The application of the 
theory to individual SNRs and their known synchrotron 
emission (jVoUd 12003 iBerezhkd [20051: IBerezhko fc Vcilk I 
|2006|) has demonstrated its capability of explaining the 
observed SNR properties and in calculating new effects 
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like the extent of magnetic field amplification which leads 
to the concentration of the highest-energy electrons in a 
very thin shell just behind the shock. The theory should 
therefore be able to explain major characteristics of the 
observed GCR spectrum up to an energy of ~ 10 15 eV 
(jBerezhko k Ksenofontovl Il999f ) under the assumption 
that the ambient magnetic field strength considerably 
exceeds typical ISM values. 

Here we apply kinetic nonlinear theory assuming a 
time-dependent amplified magnetic field, consistent with 
multi- wavelength evidence from individual objects, in or- 
der to approximately calculate the spectra of CRs pro- 
duced in Galactic SNRs. It is shown that these spectra 
are approximately consistent with the existing measure- 
ments of GCR spectra up to an energy of 10 17 eV. 

2. MODEL 

Our nonlinear model is based on a fully time-dependent 
solution of the CR transport equation together with the 
gas dynamic equations in spherical symmetry. Since all 
relevant equations, initial and boundary conditions for 
this model have already been desc ribed in detail else- 
where (e.g. IBerezhko k Volkl l2004f ) . we do not present 
them here and only briefly discuss the most important 
aspects below. 

The basic scenario is that of a point explosion, where 
the supernova explosion ejects an expanding shell of mat- 
ter with total energy E sn and mass M e j into the sur- 
rounding ISM. Due to the streaming instability CRs effi- 
ciently excite large-amplitude magnetic fl uctuations up- 
stream of the SN shock (e.g. iBelll 1 19781) . Since these 
fluctuations scatter CRs extremely strongly, we assume 
the CR diffusion coefficient to be as small as the Bohm 
limit k(p) = pvc/(3ZeB), where e and m are the proton 
charge and mass, v and p denote the particle velocity 
and momentum, Z is the particle charge number, B is 
the magnetic field strength, and c is the speed of light. 
Note that due to the dependence k cx vp/Z the shock- 
produced ion distributions depend on rigidity w ith the 
^utoff energy e max oc Z (e.g. IBerezhko et aT1ll996l ). 

If Bism is the preexisting field in the surrounding 
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ISM, then the Bohm limit implies that the instability 
growth is restricted by some nonlinear mechanism to 
the level SB ~ Bismj where SB is the wave field. An 
earlier attempt for a nonlinear descriptio n of the mag- 
netic field evolution i n a n umerical model (|Lucek k, Belli 
l2000t iBell k, Lucek 1 1200 ll ) and more recent instability 
calculations (Be 111 120041) concluded that a considerable 
amplification to what we call an effective magnetic field 
B ^> Bism should occur. The Bohm limit is then ex- 
pected to refer to this amplified field. We adopt this 
point of view here. 

From an analysis of the synchrotro n spectrum of 
SN 1006, Cassiopeia A, Tycho's SNR (see lBerezhkoll2005l 
for review) such a strong magnetic field amplification can 
only be produced as a nonlinear effect by a very effi- 
ciently accelerated nuclear CR component. The same 
large effective magnetic field is required by the compar- 
ison of our selfconsistent theory with the morphology of 
the observed X-ray synchrotron emission, in particular, 
its spatial fine structure. 

In fact, for all the thoroughly studied young SNRs, the 
ratio of magnetic field energy density Bq/8it in the up- 
stream region of the sh ock precursor to the CR pressure 
P c is about the same (jVolk et al.ll2005f) : B%/(8itP c ) « 
5 x 10~ 3 . We note that this amplified magnetic field Bo 
exceeds the typical ISM value -Bism ~ 5 mG during much 
of the early evolution. Our present calculations assume 
that this is the case for most of the evolution of the SNR 
and to this extent the results of our calculations below 
are insensitive to the concrete value of Bism except in 
the very late SNR evolutionary phase. 

The number of suprathermal protons injected into the 
acceleration process is described by a dimensionless injec- 
tion parameter rj which is a fixed fraction of the number 
of ISM particles entering the shock front. We adopt here 
a value rj ~ 10~ 4 , which is consis tent with the theoreti- 
cal expectation ()Volk et alj 120031 ) and which is close to 
the values determined individually for the three objects 
SN 1006, Tycho's SNR, and Cas A. 

The injection rate of ions heavier than protons can 
not be calculated with the required precision. Therefore 
we use ion injection rates which provide ion-to-proton 
ratios as observed in the GCRs at an energy of 1 TeV. 
The physical factors (ion injection rate and acceleration 
efficiency) which determine this ratio were discussed by 
iBerezhko fc Ksenofontovl (| 19991 ). 

The overall CR spectrum A^e^Tsn) is formed dur- 
ing the active period of SNR evolution which lasts up 
to the time Tsn when the SN shock becomes too weak 
to accelerate efficiently a new portion of freshly injected 
particles. After their release from the parent SNRs the 
accelerated CRs occupy the confinement volume more or 
less uniformly with an intensity J(efc) oc T csc (R)N(ek) 
where r esc (R) is the mean residence time, R = pc/(Ze) 
is rigidity, p is momentum, e = Amc 2 + ek and tk are 
total and kinetic particle energy, respectively. 

3. RESULTS AND DISCUSSION 

We use the values = 10 51 er g f° r the explosion 
energy and M e j = 1AM@ for the ejecta mass which are 
typical for SNe la in a uniform ISM. Note that the main 
fraction of the core collapse SNe has relatively small ini- 
tial progenitor star masses between 8 and 15 Mq which 
therefore do not significantly modify the surrounding 
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Fig. 1. — GCR intensities at the Solar system as a function 
of particle kinetic energy. Curves represent the CR spectra cal- 
culated for a proton injection rate r\ = 10~ 4 , ion-to-proton ratios 
which provide a fit for the data at = 1 Te V, and for p, = 0.75 . 
Experim ental data obtained in the CAPRICE (Boczio et al.H2003l ). 
ATIC-2 llPanov et al.ll20Qif3). JA CEE QAsakimori et alJ^OOa Tand 
KASCADE (IAntoni*e t al. 2005) experiments are shown as well. 



ISM thro ugh the main sequence wind of the progenitor 
star (e.g. Abbott 1982). SNR evolution in this case is 
similar to that of SNe la. 

The active phase of the average SNR as CR source 
was assumed to last until an age of T§n ~ 10 5 
(|Berezhko et al.ll2003D . 
In Fig.l we present the calculated intensities J(ek) X 



yr 



ef of protons (H), Helium, three groups of heavier nu- 
clei, and "All particles" as a function of kinetic energy, 
as solutions of the nonlinear equations. Here we have 
used Tosc R~^, with /i = 0.75. The results of the re- 
cent experiments CAPRICE, ATIC-2 , JACEE and KAS- 
CADE, which in our view are the most reliable ones, 
agree quite well with this theoretical calculation up to 
the energy ~ 10 17 eV. One can see that the theory fits 
the existing data in a satisfactory way up to the energy 
6^ 10 17 eV. The main exception - and difficulty - is 
the Helium spectrum as measured in the recent ATIC-2 
balloon experiment which is noticeably harder than the 
proton spectrum, in contrast to the theoretical expecta- 
tion. 

The second difficulty for the present theory are the cal- 
culated very hard overall CR source spectra which, say 
for protons, have the form N oc e^ 1 ' 9 for < 10 15 eV. 
This spectrum is noticably harder than the source spec- 
trum N oc £fc 21 , deduced in the framework of our pre- 
ferred CR propagation model which incl udes a selfcon - 
sistent halo within a Galactic wind (jPtuskin et al.lfl997l ). 
Note that this latter model explains well the existing 
data, except for the observed low GCR anisotropy which 
is at the moment a common proble m for all CR diffusio n 
models without re-acceleration (e.g. lPtuskin et aTT 2006). 
even though the anisotropy is presumably determined by 
the local structure of the Galactic magnetic field and may 
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Fig. 2. — All-particle GCR intensity as a function of total par- 
ticle energy. The dashed line represents the Galactic component, 
which is the all-particle spectrum from Fig.l. The dash-dotted 
line represents the assumed extragalactic component. Experimen- 
tal data o btained in the ATI C-2, JAC EE, KASCADE, A keno 
- AGAS A lITakeda et ap 120031) . HiRes l|Abbasi et alj 120051 1 and 
Yakutsk (Egorova et al. 2004) experiments arc shown as well. 

deviate significantly fro m the global chara cteristics of 
the propagation model (IPtuskin et al.l 119971 ) . However , 
IPtuskin fc Zirakashvilil (2003, 2005) have argued that al- 
ready in the middle Sedov phase nonlinear dissipation 
will increasingly reduce the field amplification and parti- 
cle scattering below the Bohm limit. High-energy parti- 
cles will then increasingly leave the SNR which can only 
continue to accelerate lower and lower energy particles. 
Although there is a need for a more detailed analysis of 
this effect, it should lead to some softening of the over- 
all GCR spectrum as the result of escape at high and 
continued acceleration at lower energies in older SNRs. 
In this sense we do not consider the calculated overall 
source spectra a problem for the arguments put forward 
here. We will study this question in detail in a subse- 
quent paper. 

According to Fig.l the knee in the observed all-particle 
GCR spectrum has to be attributed to the maximum en- 
ergy of protons, produced in SNRs. The steepening of 
the all particle GCR spectrum above the knee energy 
3 x 10 15 eV is a result of the progressively decreasing 
contribution of light CR nuclei with increasing energy. 
Such a scenario is confirmed by the KASCADE experi- 
ment which shows relatively sharp cutoffs of the spectra 
of various GCR spe cies at energies e max ss3Zxl0 15 eV 
(|Antoni et al.| [2005'). so that at energy e k ~ 10 17 eV the 
GCR spectrum is expected to be dominated by the con- 
tribution from the iron nuclei. 

We note that the maximum CR energy in the over- 
all CR spectrum is determined by the CRs produced 
in the very beginning of the Sedov phase. During the 
free expansion phase the SN shock produces even higher- 
energy particles as a result of the higher shock velocities 
and magnetic field strengths. However, since the mass 
in the very fast ejecta that produce them is very small, 
these particles contribute only a very steep part (tail) 
of the final overall CR spectrum, as was shown earlier 
(jBerezhko fc Volkll2004h for type la SNe. For this gen- 
eral reason the same should also be true for wind SNe. 
As a result these particles play no role for the popula- 
tion of Galact i c CRs - in contrast to the assump tion of 
iBell fc Luce"kl (|2001l) which has been followed by Hillas 
(|2006h . 



In Fig. 2 we present an all-particle spectrum which in- 
cludes two components: (i) CRs produced in SNRs and 
(ii) extragalactic CRs, protons plus 10% Hel ium, presum- 
ably produced in Active Galactic Nuclei (jAloisio et al.l 
2006). The second component has been chosen to have 
a power-law source spectrum J s oc e -2 7 above an energy 
e ~ 10 16 eV up to energies in excess of 10 20 eV. 

Compared with the source spectrum J s (e), the com- 
ponent J(e), observed in the Galaxy, is modified by two 
factors. At energies e > 10 18 eV the shape of J(e) is 
influenced by the energy losses of CRs in intergalactic 
space as a result of their interaction with the cosmic 
microwave background that leads to the formation of a 
" dip" structure at e ~ 10 19 eV, and to a GZK-cutoff for 
e > 3 x 10 19 eV (lAloisio et al.ll200(l . 

For e < 10 18 eV the spectrum J(e) is determined by 
the character of CR propagation in intergalactic space. 
Since we assume the existence of a Galactic Wind, CRs 
penetrating into the Galaxy from outside are in addition 
subject to modulation by the wind. We describe this ef- 
fect by the modulation factor / = exp(— e m /e), where the 
maximum CR energy modulated by the Galactic Wind 
is about e m = 10 17 eV dVolk fc Zirakashvilill200l . 

Cf. lAloisio et al.l |2006), presenting data of the 
AGASA, Yakutsk and HiRes detectors in Fig. 2, we shift 
the energies by a factor of A = 0.9, 0.85 and 1.2 respec- 
tively, so that the measured CR fluxes agree with each 
other. 

According to Fig. 2 the calculated GCR spectrum is in 
reasonable agreement with the existing data. It leads 
us to the following conclusion: if the observed CR spec- 
trum at energies e > 10 18 eV is indeed dominated by 
the contributio n from extragalactic sources [the so-called 
"dip-model" of lAloisio et~atl (|2006t )]. then we do not need 
any other Galactic source population except SNRs, as 
calculated above. However, if the extragalactic sources 
produce in reality a much harder spectrum J s oc e~ 2 
(the so-called " ankle-model" ) then their contribution be- 
comes dominant only at energies e > 10 19 eV. Therefore, 
to fit the observed GCR spectrum an additional Galactic 
source population is required whose contrib ution is essen - 
tial in the energy range 10 17 < e < 10 19 eV (|Hillasll2006h . 
It co uld possibly result from CR rea cceleration processes 
(e.g. lBerezhko fc Ksenofontovlll999f ). for example, in the 
interaction regions of the Galactic Wind induced by the 
spira l structure in the Galactic Disk (|Volk fc Zirakashvilil 
l2004f) . For further Galactic particle sources. see lHorandell 

pool . 

We note that these two scenarios for the extragalactic 
component predict a very different CR chemical compo- 
sition above about 3 x 10 17 eV. Since within the range 
3Z x 10 15 < e < Z x 10 17 eV reacceleration produces a 
power law tail of the CR spectrum originally produced 
in SNRs, it is clear that the observed CR spectrum is 
expected to be dominated by the iron contribution at 
energies 10 17 < e < 10 19 eV. It is very much different 
from what is expected within the dip-model. 

In order to illustrate the CR chemical composition, 
expected in the latter case, we present in Fig. 3 the 
mean logarithm of the GCR atomic number A. At en- 
ergies £k < 1 TeV A(e) increases with energy and for 
10 3 < e k < 10 6 GeV we have <lnA> « 6. Due to the 
dependence of the maximum energy of CRs produced in 
SNRs, e max cx Z, on the charge number Z and therefore 
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Fig. 3. — Mean logarithm of the CR nucleus atomic number as a 
function of particle kinetic energy. Experimental data obtained 
in the ATIC-2 , JACEE, KASCADE (QGSJET and SYBYLL ; 
IHorandel|[2005ll. as well as in the HiRes (QGSJET; IHorandelll2003T ) 
and Yakutsk (Ivanov et al. 2003) experiments are shown. 

on A, CRs become progressively heavier above the pro- 
ton cutoff energy efc ss 3 x 10 15 eV. The mean atomic 
number of the CRs produced in SNRs goes towards the 
value A = 56 as the energy approaches 10 17 eV. However, 
already at tk > 10 16 eV the contribution of extragalactic 
CRs becomes essential. Therefore < In A > reaches its 
peak value at efc ~ 2 x 10 16 eV and then diminishes with 
increasing energy to the value 0.1 x In 4. 

The calculated atomic number agrees reasonably well 
with the existing data up to a particle energy of ~ 
10 16 eV. At higher energies A(ek) goes down with energy 
similar to the Yakutsk and HiRes data. Quantitatively 
the calculated value of < In A > agrees well with the 
HiRes data. Since the existing data demonstrate that 
above 10 17 eV CRs become lighter with energy, this could 
be considered to favor the dip-scenario. 

4. SUMMARY 

Magnetic field amplification leads to a considerable in- 
crease of the maximum energy of the CRs accelerated in 
SNRs. Calculations performed within nonlinear kinetic 
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theory demonstrate that the expected GCR spectrum, 
produced in SNRs, fits the existing GCR data in a sat- 
isfactory way up to the energy 10 17 eV. The first knee 
in the observed all-particle GCR spectrum is attributed 
to the maximum energy of protons that are produced in 
SNRs. The steepening of the all-particle GCR spectrum 
above the knee energy 3 x 10 15 eV is then the result of 
a progressive depression of the contribution of light CR 
nuclei with increasing energy. Such a scenario is con- 
firmed by the KASCADE experiment which shows rel- 
atively sharp cutoff s pectra of GCR spe cies at energies 
emax ~ 3Z x 10 15 eV (jAntoni et al.|[2005l ). A difficulty of 
the present computation are the resulting very hard over- 
all source spectra from SNRs. We expect however that 
inclusion of wave dissipation and particle escape from 
SNRs at later phases of their evolution will rectify this 
extreme. 

If the Extragalactic CR component has a spectrum 
J s oc e~ 2 7 so that it dominates the observed GCR spec- 
trum already at e > 10 18 eV, then GCRs at lower energies 
are produced in SNRs without significant contribution 
from other Galactic CR sources. The GCR spectrum up 
to the energy e = 10 17 eV is dominated by the contribu- 
tion of SNRs, whereas at e > 10 18 eV GCRs are predom- 
inantly Extragalactic. If the Extragalactic source spec- 
trum is much harder, J s cx e~ 2 , then the transition from 
the Galactic to the Extragalactic component is expected 
at higher energy e ~ 10 19 eV. Since the expected GCR 
chemical composition at 10 17 — 10 18 eV is very different 
in these two cases, an experimental study of GCR com- 
position could discriminate them. The existing HiRes 
and Yakutsk data favor the first, so-called dip scenario. 
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